Statics and dynamics of atomic dark-bright solitons in the presence of delta-like impurities 
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Adopting a mean-field description for a two-component atomic Bose-Einstein condensate, we study the stat- 
ics and dynamics of dark-bright solitons in the presence of localized impurities. We use adiabatic perturbation 
theory to derive an equation of motion for the dark-bright soliton center. We show that, counter-intuitively, 
an attractive (repulsive) delta-like impurity, acting solely on the bright soliton component, induces an effective 
localized barrier (well) in the effective potential felt by the soliton; this way, dark-bright solitons are reflected 
from (transmitted through) attractive (repulsive) impurities. Our analytical results for the small-amplitude oscil- 
lations of solitons are found to be in good agreement with results obtained via a Bogoliubov-de Gennes analysis 
and direct numerical simulations. 



I. INTRODUCTION 



The physics of atomic Bose-Einstein condensates (BECs) 
01 Q] has offered the possibility of the study of purely non- 
linear phenomena in the mesoscopic scale. In particular, 
there has been a vast amount of research efforts devoted to 
the study of macroscopic nonlinear excitations of BECs (see, 
e.g., Refs. OrtZI] for reviews on this topic). In that regard, 
of particular interest are the so-called matter-wave solitons, 
of either the bright J6|] or the dark [7] type, that can be sup- 
ported in BECs with attractive or repulsive interactions, re- 
spectively. Nevertheless, these types of solitons may coex- 
ist in multi-component condensates with repulsive interac- 
tions (see, e.g., Refs. |@] and H for relevant work in two- 
component and spinor condensates, respectively): such, so- 
called dark-bright (DB) solitons exist due to the fact that the 
dark-soliton component creates, through the inter-species in- 
teraction, a trapping mechanism that allows the bright soli- 
ton component to be formed (even though this is not possible 
in repulsive BECs). Dark-bright solitons have been studied 
extensively in other contexts, such as nonlinear optics II 1 CHI 
and the theory of nonlinear waves ifTTl l. Furthermore, they 
have recently been analyzed in discrete settings |T3l . while 
higher-dimensional generalizations — namely, vortex-bright- 
soliton structures — were studied as well 11311 . Importantly, 
dark-bright solitons have been observed in experiments con- 
ducted both in optics 11141 Il5ll and, more recently, in BECs 

[EMI. 

On the other hand, the interaction of solitons with local- 
ized impurities is a quite general and fundamental problem 
that has attracted much attention in the theory of nonlinear 
waves [20] and solid state physics 02U 12211 . In this context, the 
interaction of either bright or dark solitons with <5-like impu- 
rities has been investigated in the framework of the nonlinear 
Schrodinger (NLS) equation (see, e.g., Refs. ||23T53l . while 
relevant studies have also appeared in the physics of atomic 
BECs (see, e.g., ll27l - [3lll ). In the latter context, localized im- 
purities may easily be created as sharply-focused far-detuned 
laser beams and can be used to manipulate matter-wave dy- 
namics (see, e.g., Ch. 17 of Ref. yfl), while they have already 
been used in experiments for the creation of solitons If32l l33ll . 



Nevertheless, to the best of our knowledge, the problem of the 
interaction of matter-wave dark-bright solitons with localized 
impurities has not been addressed so far. 

In this work, we aim to study this problem in the frame- 
work of mean-field theory. More specifically, we consider a 
quasi one-dimensional (ID) two-component repulsive BEC, 
composed by two hyperfine states of the same alcali species 
(as, e.g., in the experiments of Refs. 0171 - 1 1910 — a system 
that can be approximated by two coupled ID Gross-Pitaevskii 
equations (GPEs) (see, e.g., Refs. OIHlZI]). We assume that 
both components are confined by the usual harmonic trap, 
while an additional small-amplitude localized (<5-like) impu- 
rity potential is also incorporated in both components. We 
employ the hamiltonian approach of the perturbation theory of 
matter-wave solitons (see, e.g., Ref. ]7|]) to study analytically 
the adiabatic dynamics of DB solitons supported in the sys- 
tem. This way, we derive an effective equation of motion for 
the DB-soliton center. We find that if the impurity potential 
acts solely on the bright-soliton component then the soliton- 
impurity interaction is effectively repulsive (attractive) for a 
genuinely attractive (repulsive) impurity. This behavior is in a 
sharp contrast with the one corresponding to dark solitons in 
single-component condensates: there, the nature of the dark 
soliton — impurity interaction is the same as the type of the 
impurity (i.e., repulsive/attractive for repulsive/attractive im- 
purities, respectively) 12711 . 

We study the statics and dynamics of solitons near the fixed 
points of the effective potential associated to the above men- 
tioned equation of motion using both our analytical approach 
and numerical simulations. We also perform a Bogoliubov - 
de Gennes (BdG) analysis to investigate the excitation spec- 
tra of the DB-solitons proper and study their stability. Where 
appropriate, we find a very good agreement between the ana- 
lytical predictions and the numerical findings, e.g., the char- 
acteristic frequencies obtained by the equation of motion and 
the eigenfrequencies of the internal modes (also known as 
"anomalous modes" [2, 5]) associated with DB-solitons. 

The paper is structured as follows. In Section 2 we present 
the model, use perturbation theory, and derive the equation of 
motion for the soliton center. In section 3, we analyze the ef- 
fective potential and forces acting on the dark-bright solitons 
and identify the most interesting case, i.e., when the impu- 
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rity acts solely in the bright-soliton component. Section 4 is 
devoted to a systematic comparison of our analytical findings 
with simulations, including results of the BdG analysis. Fi- 
nally, in section 5 we summarize our conclusions. 



II. MODEL AND ANALYTICAL CONSH)ERATIONS 
A. Setup 

We consider a two-component elongated (along the x- 
direction) repulsive BEC, composed of two different hyper- 
fine states of the same alkali isotope. Assuming that the trap is 
highly anisotropic, with the longitudinal and transverse trap- 
ping frequencies being such that uj x -c oj±, we may describe 
this system by the following two coupled GPEs iHlj]: 



ihdtipj = ( ~2^ 9 x + Vj(x) - (ij 5j*#fc| 2 J i'P (!) 

Here, , 4>j(x,t) (j — 1,2) denote the mean-field wave func- 
tions of the two components (normalized to the numbers of 
atoms Nj = J_°° | tpj \ 2 dx), m is the atomic mass, fij are the 
chemical potentials, g.^ = 2haj±ajk are the effective ID cou- 
pling constants, ajk denote the three s-wave scattering lengths 
(note that a%2 — a 2 i) that account for collisions between 
atoms belonging to the same (cijj) or different (djk,j 7^ fc) 
species, and Vj(x) represent the external trapping potentials. 

We assume that both components are confined by the usual 
harmonic trap, namely V(x) = (l/2)muj 2 x 2 , while an addi- 
tional localized "impurity" potential, which may be created by 
a far-detuned laser beam, is also present. If such an impurity is 
strongly localized, one may theoretically approximate its spa- 
tial profile by a ^-function; thus, the trapping potentials for the 
two components can be described as: Vj(x) = V(x) + bj5(x), 
where bj are the barrier amplitudes in each component. Note 
that for a blue- or red-detuned laser beam, the impurity poten- 
tial can either repel (bj > 0) or attract (bj < 0) the atoms of 
the respective component of the condensate. 

We examine the case where the two-component BEC un- 
der consideration consists of two different hyperfine states of 
87 Rb, such as the states |1, —1) and |2, 1) used in the experi- 
ment of Ref. Q, or the states |1, -1) and |2, -2) used in the 
experiments of Refs. 1 171 - 1 1911 . In the first case the scattering 
lengths take the values an = 100. 4ao, a\ 2 = 97.66ao and 
(222 = 95.00ao, while in the second case the respective values 
are an = 100. 4ao, a\2 — 98.98ao and 022 = 98.98ao (where 
ao is the Bohr radius). In either case, the scattering lengths 
take approximately the same values, say o^- a, which is 
what we will assume hereafter. Thus, measuring the densities 
\ipj\ 2 , length, time and energy in units of 2a, a± — \/ H/uJ±, 
wj 1 and hw_\_, respectively, we may cast Eqs. (UJ into the fol- 



lowing dimensionless form, 

id t u d = - ^d 2 u d + V d {x)u d 

+ (M 2 + Kl 2 - fj)u d , (2) 

id t u b = - -dlu b + V b (x)u b 

+ (M 2 + M 2 - n - A)u b . (3) 

In the above equations, we have used the notation ipi = Ud 
and %p2 = u b , indicating that the component 1 (2) will be 
supporting a dark (bright) soliton. Notice that the respective 
normalized chemical potentials read fi% = fx^ = fi and /X2 = 
/Xb = /j, + A, and below we will assume that fi d > 
A = — |A| < 0). Finally, the external potentials in Eqs. (f2]i- 
(O take the form 

V d {x) = V{x) + hd(x) = ^n 2 x 2 + biS(x) (4) 
V b (x) = V(x) + b 2 5{x) = ^n 2 x 2 + b 2 S(x), (5) 

where ft — ui x /uj± and bi, b 2 are the normalized trap strength 
and barrier prefactor strength, respectively. Below, both of 
these parameters will be considered to be small, i.e., il ~ b <C 
1. 

Before proceeding further, it is necessary to consider at 
first the effect of the impurity on the Thomas-Fermi (TF) 
cloud carrying the dark soliton. According to the analysis 
of Ref. 12711 . the TF density near the trap center (where the 
impurity is located) can be approximated as: 

| WTF | 2 w M - 2y/Jif(x), (6) 

f(x) = i^ 2 + ^exp(-2v7iM), (7) 

where f(x) is considered to be small with respect to the chem- 
ical potential fi. The first term in the right-hand side of Eq. (jT) 
accounts for the unperturbed TF density (in the absence of the 
impurity); on the other hand, the second term actually approx- 
imates the delta-like impurity, which creates in the TF density 
a localized dip (hump) for b\ < (61 > 0); the latter has obvi- 
ously a discontinuous derivative at x — due to the matching 
conditions at x — (see details in Ref. Il27ll ). 



B. Perturbation theory 

We assume that the dark soliton is on top of a modified TF 
cloud, as described by Eqs. ©-dill. Accordingly, the density 
\ud\ 2 in Eqs. ©-Q is substituted by \u d \ 2 — > |uTF| 2 |«d| 2 - 
Furthermore, introducing the transformations t — > fit, x — > 
^JJlx, \u b \ 2 — > /i _1 |ub| 2 , we cast Eqs. ©-(O into the form: 

id t u d + \^d 2 u d - (M 2 + |it 6 | 2 - l)u d = Rd, (8) 
idtUb + ^d 2 u b - (k| 2 + k| 2 - p)ub = R b , (9) 
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where p, = 1 + A/fi, and 

R d = (2a* 2 )- 1 [2(1 - \u d \ 2 )V(x)u d + V (x)d x u d \ 



bin' 



-1/2 



(l - \u d \ 2 ) U d - y-d x U d 



-2\x\ 



Rb 



(1 - \u d \ 2 )V(x)u b + b 2 ^8{x)u b 



b U ^ 2 \u d \ 2 u b e- 2 ^ 



(10) 



(11) 



with V'(x) = dV/dx. Equations ©-(O can be viewed as 
a system of two coupled perturbed NLS equations, with per- 
turbations given by Eqs. (TTOb-dTTb. In the absence of the per- 
turbations (f2 = 0, bi t 2 = 0), and considering the boundary 
conditions \u d \ 2 — > 1 and \u b \ 2 — > as \x\ — > oo, the NLS 
Eqs. (O-© possess an exact analytical DB soliton solution of 
the following form (see, e.g., Ref. H]): 

u d (x,t) = cos^tanh [D(x — xo(t))] + isin</>, (12) 
Ub(x,t) = r/sech [D{x — xo(t))] exp [ikx + i9(t)] , (13) 

where <f> is the dark soliton's phase angle, cos <f> and rj rep- 
resent the amplitudes of the dark and bright solitons, D and 
xo(t) denote the width and the center of the DB soliton, while 
k = D tan cf) — const and 0(t) are the wavenumber and phase 
of the bright soliton, respectively. The above parameters of the 
DB-soliton are connected through the following equations: 



D 2 

x 

8{t) 



cos <p — 1] 
D tan (/>, 
1 



(D 2 - k 2 )t+{A/fi)t, 



(14) 
(15) 

(16) 



where xq is the DB soliton velocity. Notice that the amplitude 
r\ of the bright soliton, the dark-soliton component's chemi- 
cal potential fi, as well as the width D of the DB-soliton are 
connected with the number of atoms of the bright soliton by 
means of the following equation [for the variables appearing 
in Eqs. ©-©I: 



N h = 



+0 °, |2 , 
\Ub\ ax = 

oo 



D 



(17) 



Let us now assume that the DB-soliton evolves adiabatically 
in the presence of the small perturbation, and employ the 
Hamiltonian approach of the perturbation theory for matter- 
wave solitons (see, e.g., Refs. 14j,|7]]) to study the DB-soliton 
dynamics. We start by considering the Hamiltonian (total en- 
ergy) of the system of Eqs. (O-©, when the perturbations are 
absent (R d = = 0), namely, 



E = - 
2 



£dx, 



£ = \d x u d \ 2 + \d xUb \ 2 + (\u d \ 2 + \ Ub \ 2 -i) 2 

- 2(A/ A1 )K| 2 . 



(18) 



The energy of the system, when calculated for the DB-soliton 
solution of Eqs. (fT2l -([T3l. takes the following form: 



4 n /l , 9 A 

E = -D 3 + x[t:D 2 sec 2 cb-- 

5 \2 ji 



X,^. (.9) 



Since we have considered an adiabatic evolution of the DB 
soliton, we may assume that, in the presence of the per- 
turbations of Eqs. dTob- (fTTT >. the DB soliton parameters be- 
come slowly-varying unknown functions of time t (see, e.g., 
(7[]). Thus, the DB soliton parameters become (f> — > <fi(t), 
D — > D(t), and, as a result, Eqs. (fl4l-([T5]> read: 



D 2 (t) = cos 2 4>(t)-^ X D(t), 
x (t) = D(t)taxi(t>(t), 



(20) 
(21) 



where we have used Eq. H7i . The evolution of the parameters 
4>{t), D(t) and xo(t) can be found by means of the evolution 
of the DB soliton energy. In particular, employing Eq. 1191 . it 
is readily found that 



dE ■ 9 
— = ADD 
dt 



\D sec 2 4>(D + D(j> tan < 



(22) 



On the other hand, using Eqs. (HJ)-© and their complex con- 
jugates, it can be found that the evolution of the DB soliton 
energy, due to the presence of the perturbations, is given by: 



dE 

— = -2Re 
dt 



{R* d d t u d + Rld t u b ) dx\ , (23) 



where asterisk denotes complex conjugate. Substituting R d 
and R b [cf. Eqs. (TTOb-dTTbl into Eq. ( 1231 and evaluating the 
integrals, we obtain from Eqs. d2Qb . (fJTJ, ( f22l and Eq. d23b 
a system of three equations for the evolution of the soliton 
parameters cj>(t), D(t) and xo(t). This system is linearized 
around its fixed point (see details in the Appendix) and, in 
the physically relevant case of sufficiently small x, leads to 
the following equation of motion for the small-amplitude dis- 
placement Xq of the soliton position from the trap center: 



Xn = 



dVcs 
dX - 



(24) 



where we have used the variables used in Eqs. <J2J-<l3j) - The 
effective potential in Eq. (l24l is given by 



V cB (X ) = ^uj 2 osc X 2 + b S cch 2 (D a X ), 



(25) 



where the oscillation frequency u osc and the parameter b are 
respectively given by: 




3a/ 1- 
1 



(i)' 



(26) 



8A)A> + x(2A) - Do) 
x [2 (1 + 2D 2 ) h + X D bi - 3 X D 2 b 2 ] , (27) 

and Dq and Dq are constants of order 0(1) (see Appendix). 
Equation (l24l has the form of an equation of motion for a 
classical particle, with the coordinate X , moving in the effec- 
tive potential V c ff . Note that in the absence of the impurities 
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[bi = b 2 = 0, i.e., b = in Eq. (l25Vl. Eq. (|24| recovers the 
results of Ref. @]: according to this work, a DB-soliton oscil- 
lates in a harmonic trap of strength £1 with the frequency uj osc , 
given in Eq. d26b ; this frequency depends on the parameter \, 
i.e., the number of atoms Nf, of the bright soliton [see the def- 
inition of x in Eq. (1191) 1, Below, we analyze the more general 
case, studying the effect of the impurities on the statics and 
dynamics of DB-solitons. 



III. THE EFFECTIVE POTENTIAL AND FORCES 



The part of the effective potential (1251 1 induced by the im- 
purities consists of three different terms, as seen by the ex- 
pression of the constant b in Eq. d27l ). Taking into regard that 
Dq and Dq are of order 0(1) and the parameter \ is small (as 
mentioned above - see also the Appendix), it is readily ob- 
served that the sign of the parameter b is mainly determined by 
the leading-order term, oc 2b\ (l + 2_Dq), in Eq. d27l ). Thus, 
it is clear that the term oc sech 2 (Z?o^o) m the effective poten- 
tial d25l l is either a localized barrier (for b\ > 0) or a localized 
well (for bi < 0). Here we should note that although Eq. (1231 
is formally valid for small a numerical investigation of the 
more general case corresponding to values of \ of order O(l) 
reveals that the nature of the potential is correctly captured by 
the above analysis. This can be understood by the fact that, 
generally speaking, increase of \ results in a decrease of Dq 
from its maximum value (which is Dq = 1) [see Eq. ( 1A41 i in 
the Appendix] and, thus, the sign of b is always determined by 
the sign of b%. 

This result suggests that the form of the effective potential 
is not significantly modified due to the presence of the bright 
soliton component, as shown in Fig. Q] Furthermore, numer- 
ical simulations in the framework of the GP Eqs. ©-(O (not 
shown here) for the DB soliton dynamics confirm the above 
picture. We have found that stationary DB-soliton states exist 
at the fixed points of the effective potential, and if these sta- 
tionary states are displaced, they perform oscillations in the 
effective potential shown in the top panel Fig.Q]or, depending 
on their initial energy, they are either reflected or transmitted 
from the effective barrier in the bottom panel of Fig. Q] This 
behavior was already described in detail in Ref. [27], where 
the interaction of matter-wave dark solitons with localized im- 
purities was studied, hence we will not discuss it further. 

Below, we focus on a quite interesting situation occurring if 
the impurity acts solely on the bright soliton component, i.e., 
bx = and b 2 ^ 0: in this case, according to Eq. d25l l, the 
forces acting on the DB-soliton, i.e., the force exerted by the 
harmonic trap, F tI , and the localized impurity, F lmv , read: 



F 



imp 



2aD sech 2 (D X ) tanh(D A" ), 
and the constant a, which is equal to b for b\ 



(28) 
(29) 



xD 2 b 2 



2 8D Q D Q + x(2D - D ) 



0, is given by: 
(30) 
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^ 0.15 
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FIG. 1: (Color online) The effective potential Q25} in the cases 
X = 0.7 [solid (blue) line] and \ — [dashed (red) line], i.e., in 
the absence of the bright-soliton component. The top and bottom 
panels correspond to b\ — 62 = —0.15 and 61 = 62 = 0.15, respec- 
tively; the harmonic trap strength is f2 = 0.1. Insets show details of 
the effective potentials in these cases near the trap center (where the 
impurities are located). 
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FIG. 2: (Color online) Top panel: A representative illustration of the 
dependence of the two forces acting on the soliton, _Ftr and Fimp 
[solid (red) line], on the soliton center Xq; the impurity is assumed 
to be attractive, i.e., 62 < 0. The dotted (blue) and dashed-dotted 
(black) lines show — F tI for \ < \ c ; in this case, the only fixed 
point is Xq = 0. The dashed (blue) line shows — F tr for \ > Xc\ 
in this case, there exist three fixed points. Note that since the profile 
Jimp remains qualitatively the same as \ changes, for simplicity of 
illustration, it is plotted only for a single value of \. Bottom panels: 
the effective potential of Eq. l |25| > [solid (blue) line], is plotted as a 
function of Xq, for x = 0.13 < Xc = 0.145 (left) and x = 1-3 > 
Xc = 0.145 (right), and it is compared to the actual potential Vb 
acting on the bright component [cf. Eq. {5j], indicated by the dashed 
(red) line. Parameter values are fi — 1 and Q — 0.1. 
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FIG. 3: (Color online) Similar to Fig. [2] but in the case of a repulsive 
barrier, 62 > 0. In this case, the only fixed point is Xq = 0, as 
shown by the profiles of the forces (top panel). On the other hand, 
the effective potential (bottom panel) is always attractive. 



In this case, it can readily be observed that the parameter a 
is negative (positive) when b > (b < 0), as the denomi- 
nator is positive for any value of x- This result is somehow 
counter-intuitive compared to the case where both impurities 
are present: now, if the impurity is repulsive then the respec- 
tive force is attractive (F; mp < 0), while when it is attractive, 
the force is repulsive (F; mp > 0). Physically speaking, this 
behaviour can be understood by the nature of atom-impurity 
and atom-atom interactions: for example, if 62 < (attractive 
impurity) then the impurity attracts atoms in the bright soli- 
ton component; nevertheless, since inter-atomic interactions 
are repulsive, attracted atoms in Ub repel atoms in the dark 
soliton component u&. If the number of atoms in the bright 
component is sufficiently large, then the effect of the attrac- 
tive impurity on the system is to create a repulsive potential 
which overcomes the attractive effect of the trap yielding a re- 
pulsive effective potential, as shown in the bottom right panel 
of Fig. 

The form of the effective forces, F tI and Fi mp > 0, is illus- 
trated in Figs.[2]and[3] there, the dependence of the forces on 
the dimensionless parameter \ [c£ Eq. ( [T9i >l (for a fixed trap- 
ping potential frequency 51) is sketched for the cases 62 < 
and &2 > 0, respectively. In either case, as \ changes, the 
profile of Fj mp (x) remains qualitatively unchanged; for this 
reason, and for simplicity of the illustration, in Figs. [2] and [3] 
we show only one curve for i*i mp . On the other hand, changes 
of x result in a much more pronounced change on F tr , since 
\ controls its slope. As a result, changes in % mav l ea d to 
the existence of one or three fixed points associated with the 
equation of motion ( 1241 ); the fixed points, Xq, can be found as 
solutions of the following transcendental equation: 



aD sech 2 (A)X£) tanh(D X^) = 0. (31) 



The possibility of existence of one or three fixed points can 
also be understood by simplifying the equation of motion ( 1241 ): 



first, we Taylor expand V e s(Xo) in Eq. ( 1251 ) around .Xo = 
(the location of the impurity) and derive from Eq. ( 1241 ) a 
simplified equation of motion for Xq of the following form: 



Xq = 



-co 2 eS X - KX, 



K = 



o- 



(32) 
(33) 

(34) 



Equation d32l represents the normal form of the bifurcation 
arising in this system; in particular, it indicates that the fixed 
point Xq — always exists; nevertheless, depending on x 
and the sign of b, a symmetry-breaking (pitchfork) bifurca- 
tion may take place; this way, two additional fixed points can 
emerge. Below we will study the linear stability of the fixed 
points and obtain characteristic oscillation frequencies ljq of 
small-amplitude motions around them. 



A. Attractive impurity 

First, we consider the case of an attractive impurity, 62 < 0; 
in this case, one or three fixed points may exist, as shown 
by the graphical representation of the forces Ftr and Fi mp as 
functions of the DB-soliton center Xq, for different values of 
the parameter x - see top panel of Fig. [2] it is observed that 
there exists a critical value of x> namely Xo f° r which F tr 
is tangent to F; mp at X = (see dashed-dotted line in the 
figure). Then, it can easily be seen that, as long as x < Xc 
there exists only one fixed point: Xq = (see dotted line in 
the top panel of Fig.|2]). On the other hand, for values x > Xc, 
there exist three fixed points (dashed line of the top panel of 
Fig. |2j. In other words, a typical pitchfork bifurcation occurs 
at the critical value Xc- the fixed point at the origin, Xq = 0, 
loses its stability and, for \ > Xc, two new stable (off-center) 
fixed points emerge. The effective potentials corresponding 
to the cases x < Xc and x > Xc are respectively shown in 
the left and right bottom panels of Fig. [2] and illustrate the 
symmetry-breaking after the bifurcation. 

The above qualitative discussion is also supported by con- 
sidering the simplified equation of motion d32l . which can 
also provide some quantitative results for the location and sta- 
bility of the fixed points, as well as the oscillatory motion of 
the DB solitons near the fixed points. Near Xq — 0, Eq. 
can be approximated by: 



Xq « — Ld ^Xq. 



(35) 



Equation ( 1351 ) describes the motion of a DB-soliton placed 
near the trap center (where the impurity is located). As long as 
Wg ff > the soliton will perform small-amplitude oscillations 
around the center with a frequency ui e g. 

Now, if x is increased, lj 2 s (which is positive for small x) is 
decreased and, at the critical point x — Xc, the effective oscil- 
lation frequency becomes oj 2 s = 0. The critical value Xc for 
which the fixed point Xq = becomes unstable (see dashed- 
dotted line in the top panel of Fig . |2]i can be determined by uti- 
lizing Eq. d33l - recall that lu osc , a and Dq in Eq. ( T33T > depend 
on the parameter x; if x is sufficiently small (an assumption 
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consistent with our previous considerations) then the equation 
uicsiXc) — leads to the approximate result: 



2Q 2 



Xc 



Q 2 



(36) 



which is in very good agreement with our numerical findings 
(see Sec. IV). Past this critical point, a soliton placed in the 
center of the trap will eventually move away from the center 
and perform large-amplitude oscillations; in this case, w 2 ff < 
and Eq. d33l will provide the growth rate of the relevant 
instability. 

As explained above, the symmetry-breaking bifurcation re- 
sults in the emergence of two new fixed points (see bottom 
right panel of Fig. [2]), which are approximately located at 
Xq = ±oj c ff/i/K. The stability of these nontrivial fixed 
points can be studied by considering small-amplitude pertur- 
bations of Eq. (ED, of the form X (t) = X$ + S(t), and de- 
riving an equation for the small-amplitude perturbations 5(t): 



5 = —U!nS 



2 2 



aDl seen 2 (D X£) 
[3 seen 2 (D X£) - 2] . 



(37) 



(38) 



Naturally, this formula applies to an fixed Xq, including 
Xq = 0, in which case it retrieves the result of Eq. (l33l >. 



B. Repulsive impurity 

Let us now consider the case of a repulsive barrier, i.e., 
£>2 > 0. In this case, the impurity-induced force acting on 
the DB-soliton is attractive, i.e., a < 0. Thus, as illustrated in 
Fig. |3] and also observed from Eq. d32l for K < 0, the only 
solution of Eq. (fJTJ is a trivial fixed point, namely Xq = 0. 
In this case, we may follow the analysis exposed above and 
study the stability of Xq = 0, as well as the small-amplitude 
oscillations around it, by means of Eqs. d35l l and d33l l. but for 
a < 0. It is expected that, at least for sufficiently small val- 
ues of x, the fixed point should be stable and solitons located 
near the trap center will perform small-amplitude oscillations. 
Nevertheless, as will be shown in the next section, the fixed 
point undergoes an oscillatory instability past a critical value 
of x, through a different mechanism. 

Below we will compare the above analytical results with 
numerical simulations. 



IV. NUMERICAL RESULTS 

In this section, we will numerically investigate the existence 
of stationary DB-soliton solutions of Eqs. iff}-©, namely 
Ud = Ud{x) and Uf, — Ub{x), located at the fixed points Xq 
obtained before. We will show that such solutions do exist 
and will subsequently study the linear stability of these states 
by means of the BdG analysis (see, e.g., Refs. 010,01). The 
latter is performed as follows: we introduce the ansatz 

u d (x,t) = U d {x) + e [a(x)e iut + b*(x) e - iut ] , (39) 

u b (x,t) = U h (x) + e [c(x)e iut + d*(x)e~ lut ] , (40) 
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FIG. 4: (Color online) The top and bottom panels show, respectively, 
the real part (oscillation frequency) and the imaginary part (insta- 
bility growth rate) of the anomalous mode eigenfrequency uj a , as 
functions of x, m the case of the (sole) fixed point Xq = 0. Solid 
(blue) lines indicates uj a as obtained from Eq. i35\ . while dashed 
(red) lines show the numerical result obtained from the BdG anal- 
ysis. The regimes indicated by (1) and (3) correspond to the cases 
X < \c and x > Xc while the vertical dashed line (2) indicates the 
critical value x — Xc- 



into Eqs. (|2]i-([3]), and keeping terms of the order of the small 
parameter e, we will solve the eigenvalue problem for eigen- 
modes {a(x) , b(x) , c(x) , d(x)} and eigenfrequencies w — 
uj r + iuii (note that the stationary state is stable when = 0). 
This way, we will obtain the excitation spectrum of the rel- 
evant stationary states, including characteristic eigenfrequen- 
cies associated with the DB-solitons. Such an eigenfrequency 
is the one pertaining to the "anomalous mode" of the system 
(namely a mode characterized by a negative energy x norm 
product (2j,|5[]), which coincides with the oscillation frequency 
of the DB soliton moving near the center of the trap (similarly 
to the case of dark solitons in one-component BECs 113611 ). 
Following this procedure, we will be able to compare charac- 
teristic eigenfrequencies of the excitation spectrum with the 
oscillation frequencies wo derived in the framework of our an- 
alytical approximations. Remarkably, we will show that, gen- 
erally, there is a very good agreement between the two. 

In our numerical results below, we will fix the chemical po- 
tential to /i = 1, the normalized trap frequency to £1 = 0.1, 
and the impurity strength &2 = ±0.15 (for the repulsive and 
attractive cases, respectively). We should also note that, in 
the numerics, we have approximated the ^-profile of the im- 
purity potential by the function f(x) = 10 sech 2 (20a:). Other 
parameter values produced results qualitatively similar to the 
ones that will be presented below. 



A. Attractive impurity: hi < (a > 0) 

1. Fixed point at Xq = 

The analytical result of Eq. ((35), namely the dependence 
of the oscillation frequency ujq on x, is shown in Fig. |4] [see 



Re(aj) 
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• , it ■ • — 
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FIG. 5: A sketch showing the path of the anomalous mode eigenfre- 
quency oj a in the excitation spectrum, as the parameter \ ls varied, in 
the case of the fixed point X£ — 0. As \ ls increased, ui a moves to- 
wards the zero eigenfrequency of the Goldstone mode, collides with 
the latter, and an imaginary eigenfrequency pair emerges. The labels 
(1), (2) and (3) correspond to \ < Xc X = Xc and X > Xc> see a l so 

Fig.a 



solid (blue) lines]. On the other hand, in our simulations, we 
first confirmed the existence of a stationary DB-soliton state 
located at x = 0, and then determined its excitation spec- 
trum. The anomalous mode associated with the DB-soliton 
was found to have an eigenfrequency u> a , which is almost 
identical to ljq [see dashed (red) lines in Fig. 0). Figure [4] 
clearly illustrates the emergence of the pitchfork bifurcation, 
occurring at Xc = 0.145 [see vertical dotted line labeled by 
(2)]; the regimes (1) and (3) correspond to the cases x < Xc 
(one stable fixed point, Xq — 0, in the effective potential) and 
X > Xc (Xq — is unstable and two additional fixed points 
emerge). 

In order to better understand the origin of the bifurcation, 
in Fig. |5] we show the path of the anomalous mode eigenfre- 
quency u) a in the excitation spectrum. At first, i.e., for \ = 0, 
u a is located at Q/y/2, which is the approximate oscillation 
frequency of dark solitons (in the absence of the bright-soliton 
component) [0, H^]. In region (1), x i s increased and oj a 
moves towards the origin. When \ — Xc, u a collides with 
the zero eigenfrequency of the Goldstone mode - see region 
(2) in the figure. This collision gives rise to the emergence 
of an imaginary eigenfrequency pair, which characterizes the 
system as long as \ > Xc - see region (3). The picture shown 
in Fig.|5]complements the bifurcation diagrams of Fig. [4] with 
the regions (l)-(3) being in correspondence to each other; see 
also for a discussion of the relevant bifurcation phenomena in 
Hamiltonian systems, the recent exposition of ll37ll . 

We have also studied numerically the manifestation of the 
above mentioned instability of a stationary DB-soliton (ini- 
tially located at x = 0), by using this state as initial condition, 
and numerically integrating Eqs. (O-©. Note that to trigger 
the onset of the instability, a small random perturbation [of 
O(10 -3 ))] was added to the initial condition. The result is 
illustrated in Fig. [6j where the time evolution of a stationary 
DB-soliton is shown, for x = 1-35 > Xc- As seen in the fig- 
ure, the initially stationary DB-soliton is exponentially unsta- 
ble and eventually departs from its initial location, and starts 
performing oscillations. Notice that the soliton energy is suf- 
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t 

FIG. 6: (Color online) Contour plot showing the evolution of a 
DB-soliton, initially placed at x = 0, for an attractive impurity 
(b — —0.15), and for x = 1-35 (in this case, the fixed point Xq = 
is unstable). Top and bottom panels show the dark- and bright-soliton 
components, respectively. The dashed (white) line indicates the lo- 
cation of the impurity. 



ficiently large so that the soliton is always transmitted through 
the effective barrier located at the origin. It is clearly observed 
that the interaction of the soliton with the impurity results in a 
position shift: in fact, as the soliton moves from the one well 
of the effective double-well potential (see bottom right panel 
of Fig. |2]) to the other, it slows down at the impurity for a short 
time and, afterwards, it is transmitted to the other well. 



2. Fixed points at the minima of the effective double-well potential 

As in the case of Xfi = 0, we numerically confirmed the ex- 
istence of stationary DB-soliton states located at the nontrivial 
fixed points, and then determined their excitation spectra. In 
Fig. |7] we compare the result of Eq. d38l l with the one ob- 
tained in the framework of BdG analysis. An excellent agree- 
ment between the two is observed, up to a critical value of 
X, namely Xci = 0.32: in this regime, cjo of Eq. (f38T > [solid 
(blue) line in the top panel of Fig. |7l coincides with the real 
part of the anomalous mode eigenfrequency tu a [dashed (red) 
line]. Nevertheless, at x = Xci, the BdG analysis reveals that 
Lo a collides with the eigenfrequency u « Q, [the so-called 
Kohn (or dipolar) mode for 6 = 0], which characterizes the TF 
background [2]. This collision results in the emergence of an 
unstable excitation mode, characterized by a complex eigen- 
frequency quartet, the imaginary part of which are shown in 
the bottom panel of Fig. [7] In this case, a Hamiltonian-Hopf 
bifurcation takes place. This procedure can be better under- 
stood in the sketch shown in Fig. [8] as the parameter \ is 
increased, the anomalous mode eigenfrequency uj a is also in- 
creased, i.e., it moves to the opposite direction as compared 
to the situation shown in Fig. [5] This way, u) a eventually col- 
lides with the eigenfrequency oj sw Q , and gives rise to the 
emergence of a quartet of complex eigenfrequencies. 

The above analysis suggests that for values \ < Xci, a 
DB-soliton initially located at any of the two nontrivial fixed 
points, when displaced, will perform small-amplitude oscil- 
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FIG. 7: (Color online) Same as Fig.[4] but in the case of the nontrivial 
fixed points Xg. Solid (blue) lines indicates ui a as obtained from 
Eq. d 3 8I > , dashed (red) lines show the numerical result obtained from 
the BdG analysis, while dotted (green) line in the top panel of the 
figure indicates the eigenfrequency of the approximate Kohn mode. 
The regimes indicated by (1), (2) and (3) correspond to the cases 
X < Xd, X = Xd and x > Xa- 



enough kinetic energy to be transmitted through the effective 
barrier. This way, it moves over to the other well of the effec- 
tive double-well potential and, afterwards, the above process 
is repeated. 



B. Repulsive impurity: b 2 > (a < 0) 

In the case of a repulsive barrier impurity, we will compare 
the relevant analytical [see Eq. d35l > for a < 0] and numerical 
results (obtained by the BdG analysis). First we mention that, 
as seen in the top panel of Fig.[l0j the oscillation frequency ojq 
of the DB-soliton almost coincides with the anomalous mode 
eigenfrequency uj a , only for sufficiently small values of pa- 
rameter x- In f ac t, there exists a critical value of x, namely 
Xc2 = 0.05, where a bifurcation - similar to the one shown in 
Fig. [8]- takes place. This bifurcation results in the emergence 
of an unstable eigenmode, characterized by a quartet of com- 
plex eigenfrequencies, the imaginary part of which is shown 
in the bottom panel of Fig. [10] 
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FIG. 8: Similar to Fig. [5] but in the case of the nontrivial fixed 
points. As the parameter \ is increased, the anomalous mode eigen- 
frequency uj a moves towards the Kohn mode eigenfrequency (lo- 
cated at ui — Q) and, after the collision, a complex eigenfrequency 
quartet emerges. The regimes indicated by (1), (2) and (3) corre- 
spond to the cases x < Xci, X = Xci and \ > Xa', see also Fig. [7] 



lations at one well of the effective double-well potential. A 
direct numerical integration of Eqs. ©-(O, with initial condi- 
tion such a stationary DB-soliton state (perturbed by random 
noise), shows that this is the case indeed: a prototypical exam- 
ple is shown in two top panels of Fig. [9] where the dynamics 
of such a state is illustrated, for \ — 0.25 < \ci = 0.32. 
It is clearly observed that the DB-soliton oscillates around 
the center of one of the wells, with an oscillation frequency 
uj a Ri 0.05; this value deviates approximately 5% from the 
analytically predicted value [cf. Eq. d38tl. On the other hand, 
it is interesting to numerically investigate the manifestation of 
the predicted instability of a stationary DB-soliton state for 
X > Xci ■ Such a case, is illustrated in the two bottom panels 
of Fig. [9] where the evolution of such a DB-soliton is shown, 
for x = 0.62. It is observed that the initially quiescent DB- 
soliton starts performing small-amplitude oscillations around 
the center of one of the wells but, after a short time, it gains 




FIG. 9: (Color online) The two top panels present contour plots 
showing the evolution of a DB-soliton, initially placed at the fixed 
point x = Xg — 0.6, for a value of \ — 0.23 < %& = 0.32; in 
this case, the fixed point is stable. The two bottom panels are similar 
to the two top ones, but for the fixed point x = Xq — 1.3, for a 
value of x = 0.62 > xa', in this case, the fixed point is oscillatorily 
unstable. First- and third- (second- and fourth-) row panels show the 
dark- (bright-) soliton components. The dashed (white) line indicates 
the location of the impurity. 
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FIG. 10: (Color online) Similar to Fig. [7] but for the case of a repul- 
sive impurity (62 = 0.15) and for a DB-soliton located at x = 0. 
The critical value of parameter \ is \c2 = 0.05. 



As before, it is relevant to numerically study the manifes- 
tation of the instability in the case of a DB-soliton initially 
placed at x = 0, for \ > Xc2- A pertinent example is illus- 
trated in Fig.Qj] where the evolution of such a state is shown 
for x = 0.15. The soli ton falls into an instability, and eventu- 
ally starts to oscillate around the center of the trap. 



V. CONCLUSIONS 

We have used mean-field theory to study the statics and dy- 
namics of atomic dark-bright solitons in the presence of local- 
ized (delta-like) impurities. Our model considered a system 
of two coupled Gross-Pitaevskii equations, describing a two- 
component Bose-Einstein condensate, confined in an external 
potential composed of a harmonic trap and a pair of localized 
impurities acting on each component. 

We have employed the adiabatic perturbation theory for 
solitons to derive an equation of motion for the dark-bright 
soliton center. Our analytical approximation revealed that 




if both impurity potentials are repulsive (attractive) then the 
effective potential felt by the soliton is either a double-well 
or a harmonic potential with a localized well located in the 
trap center. Investigating the forces acting on the soliton, we 
have identified an interesting situation, which was then an- 
alyzed in detail: if the impurity potential acts solely on the 
bright-soliton component, then the impurity-induced part of 
the effective potential is either a localized barrier (for attrac- 
tive impurity) or a localized well (for repulsive impurity). This 
behavior is in sharp contrast with the one corresponding to 
the case of single-component condensates, where the nature 
(repulsive or attractive) of the soliton-impurity interaction is 
identical to the type of the impurity (repulsive or attractive) 
101. 



FIG. 11: (Color online) Similar to Fig. [6] but for the case of a repul- 
sive impurity {b-z = 0.15). 



Our numerical simulations have confirmed that stationary 
dark-bright solitons do exist at the fixed points of the effec- 
tive potential. The stability of these fixed points was studied, 
and the frequency of small-amplitude oscillations in a stable 
configuration was found analytically. We have performed a 
Bogoliubov-de Gennes analysis to study the stability of sta- 
tionary states and find their excitation spectrum. The eigen- 
frequencies of the anomalous modes were found to be almost 
identical to the analytically obtained soliton oscillation fre- 
quencies, at least for sufficiently small number of atoms of 
the bright component. In the case of unstable fixed points, the 
bifurcations (pitchfork or Hamiltonian Hopf) that give rise to 
the destabilization are identified and the growth rate of the 
perturbations are theoretically identified and corroborated by 
numerical linear stability analysis. 

An interesting direction for future studies would be a sys- 
tematic study of the scattering of dark-bright solitons from 
localized impurities of arbitrary amplitude (in the lines of the 

work in Refs. mm USE!). Furthermore, it would be in- 
teresting to study similar problems but for impurities that have 
spatial scales larger than the ones of the soliton, and investi- 
gate possible changes in the stability and dynamics. Addition- 
ally, it would be quite relevant to extend the present analysis 
(and its pertinent generalizations as per the previous points) 
to multi-dimensional settings, and study the statics and dy- 
namics of vortices in the presence of localized impurities (see, 
e.g., a relevant study but for a single-component condensate in 
Ref. HU). Such studies are in progress and pertinent results 
will be presented elsewhere. 
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Appendix A: Equation of motion for the soliton center 

Substituting R d and R b [cf. Eqs. (TTOt-dTTbl into Eq 
and evaluating the integrals, we obtain from Eq. (1231 the fol 
lowing result: 



dE 
~dt 



pT 2 sin(20) (cos 2 <t> - rf) V'{x) 

ifeiL>sin(2(/))(cos 2 cj> + 2D 2 )h 
b2X(J-~ 2 D 3 tan^sech 2 (Z?xo) tanh(Z?xo) 
\bixn~ 2 D 2 tan0(/i cos 2 4> - J 2 ), 



of atoms [ 16-19], we may approximate the above integrals as 
l x « I 2 w (2/3) sech 2 (Dxo) tanh(Dxo). This way, we ac- 
cordingly simplify Eq. (IA1I) . which together with Eqs. ( 1201 ). 
(fJTJ [and Eq. d22l l constitute a system of three ordinary dif- 
ferential equations for the unknown soliton parameters 4>(t), 
xo(t) and D(t). This system can be solved approximately 
upon linearizing around the fixed point: 



k = 0, 4 0) =0, D Q = 



1 + ' 7 



(Al) 



using the ansatz iro = Xo, <fi — <t>i an d -D 
thus obtain the following results: 



_ X 
4' 



(A4) 
Di. We 



where we have Taylor expanded the potential V(x) around the 
soliton center xo and assumed that the DB-soliton is moving 
in the vicinity of the trap center (where the impurity is lo- 
cated), xq w 0; this way, we actually deal with nearly station- 
ary DB-solitons, characterized by slow velocities, such that 
the phase angle is <f> w 0. Furthermore, I\ and I 2 in Eq. ( 1AU 
are the following integrals: 



h = 
h = 



+ 00 



— oo 



sech 4 [L>(x-a;o)]e" 2|a:| 
sech 2 [D(a;-xo)]e- 2 l a; l 



dx, (A2) 
dx, (A3) 



which can be evaluated by means of the hypergeometric func- 
tions B35I1 . Nevertheless, in the physically relevant case of suf- 
ficiently small x [cf. Eq. (U9H, i.e., when the number of atoms 
of the bright soliton is only a small fraction of the total number 



where 



= -D 



0<Pl, 



Do=[2Do + ^ 



-2 + X Do - D 

2 > 



2D 2 ) - X Do [ b 2 Do 



x sech 2 (£>oXo)tanh(D Xo), 



V = -Do 



8 2 D + x{2D - Do) 



bi 
3 



(A5) 



(A6) 
(A7) 



(A8) 



To this end, differentiating Eq. ( IA7b with respect to time once, 
and using Eq. ( IA6b , after some straightforward algebraic ma- 
nipulations, we obtain the equation of motion ( 124b for the DB 
soliton center. 
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